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ABSTRACT 

Blandford-Znajek (BZ) mechanism describes a process extracting rotation energy from a spinning 
black hole (BH) via magnetic field lines penetrating the event horizon of central BH. In this paper, 
we present a perturbation approach to study force-free jets launched by the BZ mechanism, and 
its two immediate applications: (1) we present a high-order split monopole perturbation solution to 
the BZ mechanism, which accurately pins down the energy extraction rate E and well describes the 
structure of BH magnetosphere for all range of BH spins (0 < a < 1); (2) the approach yields an 
exact constraint for the monopole field configuration in the Kerr spacetime, I = f2(l — A^), where 
is the <f>— component of the vector potential of electromagnetic field, fi is the angular velocity of 
magnetic field lines and / is the poloidal electric current. The constraint is of particular importance 
to benchmark the accuracy of numerical simulations. 

Subject headings: gravitation - magnetic fields - magnetohydrodynamics 


1. INTRODUCTION 

The Blandford-Znajek mechanism (Blandford & Zna- 
jek 1977) is believed to be an efficient way to extract 
rotation energy out of a spinning BH, which operates in 
BH systems on all mass scales, from stellar-mass BHs of 
gamma ray bursts to supermassive BHs of active galactic 
nuclei. In the past decade, general relativistic magneto¬ 
hydrodynamic simulations have shown a physical picture 
of jets launched by the BZ mechanism (Komissarov 2001, 
2004a,b, 2005; Komissarov & McKinney 2007; McKin¬ 
ney & Gammie 2004; McKinney 2005; Tchekhovskoy 
et al. 2008, 2010, 2011; Tchekhovskoy & McKinney 2012; 
Penna et al. 2013; McKinney et al. 2013), which shows 
that purely poloidal magnetic fields, Bp, appear around 
a non-rotating BH, while a spinning BH induces poloidal 
electric field Ep and toroidal magnetic field Bt, thus an 
outward Poynting flux Ep x Bp is generated, and rota¬ 
tion energy of the spinning BH is extracted in the form 
of Poynting flux. 

Independent of time-dependent simulations, Con- 
topoulos et al. (2013); Nathanail & Contopoulos 
(2014) numerically solved the general relativistic Grad- 
Shafranov (GS) equation which governs the structure of 
the stationary force-free magnetosphere around a Kerr 
BH. It is known that the GS equation is a second-order 
differential equation of A^ depending on the magnetic 
field angular velocity H and the poloidal electric cur¬ 
rent /. The GS equation was shown to be an eigen¬ 
value problem with eigenfunction A$ and eigenvalues 
I and f! to be determined. To solve the GS equation 
self-consistently, eigenvalues I and H were adjusted to 
make magnetic fields smoothly cross two light surfaces. 
In this way, they confirmed the results of previous sim¬ 
ulations for asymptotically monopole and paraboloidal 
magnetic field. But for the asymptotically uniform mag¬ 
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netic field, there is only one light surface, which is not 
sufficient for determining the two eigenvalues / and f2, 
so they argued that there are infinitely many solutions 
for the asymptotically uniform field (Contopoulos et al. 
2013; Nathanail & Contopoulos 2014; Yang et al. 2015). 
While time-dependent simulations conducted by differ¬ 
ent groups seem to converge to similar solutions (e.g., 
Komissarov (2005); Komissarov & McKinney (2007); 
Alic et al. (2012)). But different simulations did not 
precisely match, especially in the range of extreme spins, 
and there is no way to tell which simulation is more ro¬ 
bust. 

To settle down all these uncertainties and controver¬ 
sies, independent analytic works are of great value (Tan- 
abe & Nagataki 2008; Pan & Yu 2014, 2015). Though 
Contopoulos et al. (2013) criticized that the analytic per¬ 
turbation approach did not appreciate the critical role of 
light surfaces, it is no coincidence that previous perturba¬ 
tion solutions are consistent with simulations, especially 
for non-extremal Kerr BHs (e.g. McKinney & Gammie 
(2004); Tchekhovskoy et al. (2010)). In fact, perturba¬ 
tion solutions obtained are analytical and smooth, which 
certainly satisfy the requirement of magnetic field lines 
smoothly crossing the light surfaces. So perturbation 
approach is consistent with the eigenvalue approach pro¬ 
posed by Contopoulos et al. (2013), and the only dif¬ 
ference is that the eigenvalue approach implements the 
smoothly crossing requirement explicitly, while the per¬ 
turbation approach implements the requirement implic¬ 
itly. 

In this paper, we propose a perturbation approach to 
investigate the BZ mechanism analytically, which in prin¬ 
ciple enables us to obtain perturbation solutions accurate 
to any order of the BH spin. Following this approach, we 
present a high-order monopole perturbation solution to 
the BZ mechanism, and derive an exact constraint rela¬ 
tion for monopole magnetic field, which provides a cri¬ 
teria for testing the accuracy of numerical simulations. 
The paper is organized as follows. Basic equations gov¬ 
erning stationary axisymmetric force-free magnetosphere 
in the Kerr spacetime are summarized in Section 2. In 
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Section 3, we present our perturbation method and ap¬ 
ply it to the monopole magnetic field. We present results 
in Section 4 and discussion in Section 5. We conclude in 
Section 6. 

2. BASIC EQUATIONS 

Stationary axisymmetric force-free electromagnetic 
fields in the Kerr spacetime are determined by three func¬ 
tions I(Aj,). Non-trivial components of Fara¬ 

day tensor F = d^A^ — 3 U A^ in the Kerr-Schild coor¬ 
dinate (McKinney & Gammie 2004; Kerr 1963) could be 
expressed as follows (Blandforcl & Znajek 1977; McKin¬ 
ney & Gammie 2004; Pan & Yu 2014) 

Rr<j) — B(pr A(fi, r , Fg^ — F^g — ^ : (1) 

Ftr = Brt = VlA^^ , F t g = Fg t = flA^ g , (2) 

F r g = Fg r = y/^gB* , (3) 

where 

B 4> __ JE + ( 2nr - a ) sin 6A <t>,e (4 y 

AS sin 2 9 ’ 1 ’ 

and S = r 2 + a 2 cos 2 9, A = r 2 — 2r + a 2 , y/—g = S sin$. 
The energy conservation equation reads 

- n [(v/=; gF tr ), r + (y/^-gF te ) fi \ + F r gI'(A^ 

+ [(V^F^), r + (V=9F* e ),g\ = 0 , (5) 

which is also known as the Grad-Shafranov (GS) equa¬ 
tion. 

In the Schwarzschild spacetime, I = fl = 0 , the GS 
equation is simplified as 

CA+ = 0, (6) 

where the operator 

1 d / 2\ 3 l 3 l 3 

sin 9 dr \ r J dr ^ r 2 39 sin 9 39 

The Green’s function G(r, 9 ; ro, 9 0 ) to the operator C de¬ 
fined by 

CG(r, 9 ; r 0) 0 o ) = S(r - r 0 )S(9 - 0 O ), (8) 

is available (Blandford & Znajek 1977; Petterson 1974). 
Note that the solution satisfies the boundary conditions 
that G(r,9;r 0 ,9 0 ) is finite at r = 2 and approaches to 
zero at infinity. 

3. PERTURBATION METHOD 

It is notoriously difficult to exactly solve the GS equa¬ 
tion (5) due to its nonlinearity. We adopt the pertur¬ 
bation technique to attack this problem. Blandford & 
Znajek (1977) first put forward the monopole perturba¬ 
tion solution up to order of 0(a 2 ) and Pan & Yu (2015) 
extended the solution to 0(a 4 ). In this paper, we adopt 
a generalized perturbation approach, which, in principle, 
enables us to obtain perturbation solutions to any high 
order. We find that truncation at the eighth-order is 
good enough to precisely match current state-of-the-art 
numerical simulations. This method can be applied to 
other type of field configuration as well, though tedious 
calculations are inevitable. For simplicity, we focus on 
the split-monopole field throughout this paper. 


For monopole magnetic field, the Schwarzschild metric 
solution to the GS equation (5) writes as 

= 0, /o = 0, A<p = A 0 = - cos 9. (9) 

For corresponding Kerr metric solution, we define lo = 
0(71^)1^00, i = I(A,j > )\ r - ¥00 , and expand them in series, 

A<p = Aq + ci 2 A2 + u 4 A 4 + O(n^), 

w = au>i + a 3 ui 3 + a 5 W 5 + 0(a 7 ), 

* = afi + a 3 *3 + a 5 t5 + 0(a 7 ). (10) 

As done by Pan & Yu (2015), O and I could be expressed 
in terms of lo and i respectively. With the above nota¬ 
tions, the GS equation (5) could be decomposed as a set 
of linear equations 

CA n (r,0) = S n (r,0-,i n -i,u} n -i) (n = 2 , 4 , 6 ,...). (11) 

These equations could be solved one by one as follows. 
Accurate to 0(a 2 ), GS equation (5) is written as 

CA 2 (r,6) = S 2 (r,0; ii,u;i). (12) 

Before solving the differential equation, i\ and ui\ need 
to be specified by two constraints. The Znajek regularity 
condition (Znajek 1977) which requires B* to be finite 
on horizon gives 

h = sin 2 9 Q - uq^ , (13) 

and the convergence constraint that perturbation solu¬ 
tion A 2 (r, 9) should be convergent from horizon to infin¬ 
ity gives (see Pan & Yu (2015) for details) 

i\ = lo 1 sin 2 9. (14) 

Consequently, 

h = <*q sin 2 9, uq = -. (15) 
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With ii and uq specified, the source function S 2 (r, 9) is 
determined, so we can write the solution A 2 (r, 9) as an 
integral of the Green’s function 

/*00 f‘TT 

A 2 (r,9) = / dr 0 / d0 o S 2 (r o> 9 0 )G(r, 9] r 0 ,9 0 ). (16) 
J 2 J 0 

In this way, we work out all the three functions, 
A t / ) ,I(A ( j ) ) and ^1{A^) up to 0{a 2 ). Refer to Blandford 
& Znajek (1977) for explicit form of A 2 (r,9), where on 
inner/outer boundary 

A 2 (2,9) = Rsin 2 9cos9, lim A 2 (r, 9) = 0, (17) 

r—> 00 

with R = ( 67 r 2 — 49)/72. 

We could obtain solutions up to 0(a 4 ) in a similar 
manner (Pan & Yu 2015). The Znajek horizon condition 
and the convergence constraint give 

* 3 = w 3 sin 2 6 », lo 3 = 1 ^ sin 2 d+^, (18) 

and the integral of Green’s function gives 

A 4 (2, 9) = psin 2 9c\ 2 \cos9) + q sin 2 9C^ 2 \cos6), 
lim A 4 (r, 9) = 0. (19) 

r—>00 






where 
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c[ 2 ^ (cos 0) = 3 cos 0, 

/ 3 •. 1 Cj Q CC 

Cg (cos 9) = — — cos0 + — cos 3 6, (20) 


and 


9(j | 13<(3) 
_ 14 + 3920 ‘ 
~ 1.7 x 10" 2 , 

6 cr 321C(3) 
+ 9800 
~ 9.0 x 10" 4 , 


7T 2 (19385 + 774tt 2 ) | 17929399 
181440 + 7620480 

7T 2 (48955 + 918tt 2 ) 2012505017 

201600 + 508032000 

( 21 ) 


with 


(7 = 3? 


i r 


2U 2 (|) log I 


(1-f) 


^-dr = 1.3529... 


e* - 1 


-dt = 1.2020... 


( 22 ) 


where 3? designates the real part and Li 2 (x) is the second- 
order Jonquiere’s function. 

Take one step further, we could obtain quantities of 

0(a 5 ) 

* 5 = w 5 sin 2 9, oj 5 = — (1 — 0.037 sin 2 9 + 0.155 sin 4 6), 
64 

(23) 

using the Znajek horizon condition and the convergence 
constraint. 


4. RESULTS 

The energy extraction rate defined by E = 
—2n y/^gT r t d9 = 27r/ I{A^)Q.{A^)dA ( j > , where T^ v is 
energy-momentum tensor of electromagnetic field, is cal¬ 
culated accurate to 0(a 6 ) as 

E = 2tt j *(Ao)w(Ao)dAo 

= —(a 2 + 0.587a 4 + 0.355a 6 ). (24) 

24 

Note that the 6 th-order correction to the energy extrac¬ 
tion rate only depends on 11 , 3,5 and w i, 3,5 which are de¬ 
termined by Aq, 2 , 4 - See Fig. 1 for the energy extraction 
rate accurate to different orders. For comparison, cor¬ 
responding simulation results from Tchekhovskoy et al. 
(2010) are also plotted. We find that the analytic energy 
extraction rate of 0(a 6 ) better matches simulation re¬ 
sults compared to the original second-order BZ solution, 
but still underestimates energy extraction rate by ~ 30% 
for extreme spins, which indicates that even higher or¬ 
ders should be included. To improve the performance of 
analytic solutions, one way is to expand them to higher 
order of spin. We follow an alternative approach in this 
paper. We aims to find perturbation expansion series of 
faster convergence. 

Simulations and analytic works agree that, w = S2 h/ 2 
is a good approximation, where fin = a/(2+2\/l — a 2 ) is 
angular velocity of a Kerr BIT We expect a faster solution 
convergence in terms of fin instead of a. Accurate to 


a 





Fig. 1. The energy extraction rate E accurate to different or¬ 
ders. 



Fig. 2.— Distribution of angular velocity of magnetic field lines 
uj accurate to 

O(fl^), we have 

ui = -S2h + 0.43fy sin 2 0 + flj) (—5.75 sin 2 0 + 2.48 sin 4 0), 

2 (25) 

which is shown in Fig. 2. We see w ~ fln/2 for a < 0.9 
and ui decrease abruptly for a > 0.9. This turn-over was 
seen in many numerical simulations and is analytically 
confirmed here for the hrst time. In addition, take the 
extreme spin case a = 0.999 for example, we see the sub¬ 
structure that uj decreases quickly then increases slowly 
with increasing angle 0 in the range of 0 < 0 < ir/2. The 
substructure was also confirmed by previous simulations 
(Tchekhovskoy et al. 2010). 

With u expressed in terms of fin 7 the energy extraction 
rate is expressed in terms of fin as 

• 27T 

E = —(^H + F38f4 - 11.09fl|). (26) 

o 

There is a smart way to extend E to even higher order 
making use of the relation 


dE 


dE da 
da cAIh 


= 0 , 


(27) 
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where we have used the empirical constraint that dE/da 
is finite and the relation (da/dflu) a —i =0. As a result, 
eighth-order energy extraction rate E can be written as 
• 27T 

E = — (G| + 1.38G^ - 11.09G& + 6.26G|). (28) 

O 

Note that the eighth-order correction here is in principle 
different from the true correction that would be obtained 
in the usual way by solving the GS equation. The en¬ 
ergy extraction rate E accurate to different orders are 
shown in Fig. 1. It turns out that our analytic solution 
up to O(G^) perfectly reproduces the result of simula¬ 
tions for the whole range of spins. Besides, the com¬ 
monly used lowest-order approximation E = (27r/3)f2^ 
overestimates the energy extraction rate by ~ 30% for 
the extreme spin, because, it does not account for the 
abrupt decrease of w from Gr/ 2 for a > 0.9. 



5. DISCUSSION 

We find that i n = oj n sin 2 9 for n < 5 for split-monopole 
magnetic field configuration. Generally, we could prove 
the relation hold for all positive integer n. Take the n = 8 
case for example, 

Sg =sin 2 0wi(sin 2 9uir),e + sin 2 0w 3 (sin 2 9u> 5 )^ 

+ sin 2 0w 5 (sin 2 9u}^)fi + sin 2 dw 7 (sin 2 9uj\) } g 
— (* 1 * 7,0 + * 3 * 5,0 + * 5 * 3,0 + * 7 * 1 , 0 ) + 0(l/r), (29) 

where we have explicitely listed all source terms of 0 ( 1 ). 
According to the convergence requirement, the summa¬ 
tion of these 0 ( 1 ) terms vanishes, so we have = 
W 7 sin 2 9. Similarly, we find i n = ui n sin 2 9 for all posi¬ 
tive n. As a result, we have / = Gsin 2 ^Ir-xxj- Note that 
both I and £2 are functions of A^, so I/G is also a func¬ 
tion A^. We denote it as //G = f(A^). This function 
can be readily determined at infinity. It is known that 
I/GIr-j-oo = sin 2 9 and A^= — cos#, so we have 
f(—cos9) = sin 2 9 or f(x) = 1 — a: 2 . Finally, we arrive 
at an exact constraint relation 

I = n(l -A 2 ), (30) 

which holds at all radii. Recently, Penna (2015) shows 
that our constraint relation (30) is equivalent to the out¬ 
going boundary condition at infinity. 

According to this constraint relation, we show that the 
distortions of poloidal lines and variations of angular ve¬ 
locity of magnetic field lines are in fact results of force 
balance of magnetic field. We start with the explanation 
of A^q variation with BH spin: BH rotation induces a 
toroidal component of magnetic field Bx, whose hoop 
stress exerts a polar-directed force on poloidal field lines. 
In response, poloidal field lines move towards the polar 
direction (see Fig. 3), thus the fraction of magnetic flux 
increases in the polar region and decreases in the equa¬ 
torial region. In other word, jA^g | increases in the polar 
region and decreases in the equatorial region due to hoop 
stress of toroidal field generated by BH rotation. Comb¬ 
ing the constraint relation (30) with the Znajek regu¬ 
larity condition, G on horizon r = r+ can be expressed 
as 

Q ( sin 9A$fi \ / 2r + \ 

Gr — G ^ 1-AJ J \r 2 + a 2 cos 2 9) 


Fig. 3. — Illustration of distortions of poloidal field lines in¬ 
duced by BH rotation, where solid lines are poloidal field lines 
in Kerr spacetime and dashed lines are corresponding field lines in 
Schwarzschild spacetime. 


so G is determined by distortions of poloidal field lines 
and the size of BH horizon. For simplicity, we consider 
G in two regions: the equatorial region # ~ 7r/2 and the 
polar region 9 ~ 0. For the former, 


G 

G h — G 


2A^q\ 
r + ) e~ tt /2 


(32) 


The right-hand side of the above equation is proportional 
to Br/Bx, where the toroidal field Bx ~ 1/r and radial 
field Br ~ A^^/t 2 . Decreasing A^g will suppress Br. 
To keep the force balance, G must decrease accordingly 
to suppress hoop stress. Decreasing r + will increase more 
Br than Bx- To keep the force balance, G must increase 
accordingly. In real situation, both r+ and A^ q decrease 
with increasing spin a, and G/Gr is determined by the 
competition between the two. For the polar region, we 
expect vanishing distortions of poloidal field lines, due to 
vanishing toroidal field, then we have 


G|0 =o = G h /2, (33) 


for any BH spin. Similar arguments also apply to other 
field configurations, e.g., magnetic held lines moving to¬ 
wards polar region and the turning over of G/Gr were 
indeed seen in simulations of paraboloidal magnetic held 
(Tchekhovskoy et al. 2010). 


6. CONCLUSIONS 

A stationary axisymmetric force-free jet launched by 
the BZ mechanism is completely determined by three 
function, A^, I(A^) and G (A^), which are related by 
the GS equation (5). We propose an analytic approach 
to self-consistently determine the three functions. As 
a specihc application, we present a high-order split 
monopole perturbation solution to the BZ mechanism. 
The solution accurately describes the structures of 
BH magnetosphere and the energy extraction rates 
of BZ mechanism for the whole range of BH spins. 
The analytic approach also yields an exact constraint 
relation, I = G(1 — A^), for the split monopole magnetic 
held in the Kerr spacetime. The constraint relation 
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provides a useful tool for testing accuracy of numerical 
simulations. Based on this constraint relation, we find 
that the distortions of magnetic field lines and the 
variations of angular velocity D of magnetic field lines 
are results of force balance between the poloidal and 
toroidal field. 
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